In measurements on systems, an exponential function may arise whenever a quantity grows or decays at a rate proportional to its current value. Properties of many natural phenomena can be described by exponential functions. For a system that is a mixture of such decay rates, a function for the system can be expressed as multi-exponential decay function:E=ΣJ=1IA1exp(−t/T1)  (1)An inverse Laplace transform applied to the measured data having the multi-exponential behavior can be used to provide the distribution of A. Such application is not limited to one-dimensional (1D) analysis, but can be applied to two-dimensional (2D) and three-dimensional (3D) analysis. Problems that can be solved with the inverse Laplace transform method include inversion of echo trains to obtain 1D relaxation time and a 2D diffusion-relaxation distribution. There are other types of problems that can be formed in a similar single or multiple dimensional inverse Laplace transform problem. The inverse Laplace transform is applicable to other problems in which measurement data can be expressed by a multi-exponential function.
Deriving petrophysical information from nuclear magnetic resonance (NMR) logs often starts with an echo train inversion using a multi-exponential model. An NMR inversion technique involves fitting parameters and variables to the measured echo trains through an optimization procedure in the discrete space. The optimization procedure includes a series of time-consuming matrix multiplications for each measurement data sequence (echo trains), which can cause intolerable processing inefficiency in processing NMR multi-dimensional parameter space inversions.
An echo train in NMR measurement can be expressed by,
                                          M            ⁡                          (                              t                ,                Tw                ,                G                            )                                =                                    ∑                              i                =                1                            I                        ⁢                                          ∑                                  j                  =                  1                                J                            ⁢                                                                    m                    ij                                    (                                      1                    -                                          e                                              -                                                  Tw                                                      RT                                                                                          2                                ⁢                                i                                                            ⁢                                                                                                                                                                                                                                                        )                                ⁢                                  e                                                            -                                              [                                                                              1                                                          T                                                              2                                ⁢                                i                                                                                                              +                                                                                                                    (                                                                                                      D                                    j                                                                    ·                                  G                                  ·                                  γ                                  ·                                  TE                                                                )                                                            2                                                        12                                                                          ]                                                              ⁢                    t                                                                                      ,                            (        2        )            where M is the echo amplitude, mij the unknown partial porosities on a grid with index (i, j), Dj the fluid diffusivity, Tw the wait time, γ the gyromagnetic ratio, G the magnetic field gradient, TE the time spacing between echoes in an echo train, and R=T1/T2. T1 and T2 are the longitudinal and transversal relaxation times, respectively. Determination of the unknown partial porosity mij according to the measured echo train amplitude is an inversion problem. Through the process of inversion, the spin-echo decay data can be converted to a T2 distribution. This distribution represents a “most likely” distribution of T2 values that produce the echo train. With proper calibration and account for hydrogen index of the fluids in the pore space, the area under a T2 distribution curve is equal to the porosity. This distribution can be correlated with a pore size distribution when the rock is 100% water saturated. However, if hydrocarbons are present, the T2 distribution will be altered depending on the hydrocarbon type, viscosity, and saturation.
The inversion problem to determine the unknown partial porosity mij can be mathematically described by,
                              min          ⁢                      {                                                            ∑                                      p                    =                    1                                    P                                ⁢                                  (                                                            w                      p                                        ⁢                                                                  ∑                                                  q                          =                          1                                                                          Q                          p                                                                    ⁢                                                                                          ⁢                                                                                          ⁢                                                                        (                                                                                                                    ∑                                                                  i                                  =                                  1                                                                I                                                            ⁢                                                                                                ∑                                                                      j                                    =                                    1                                                                    J                                                                ⁢                                                                                                                                            m                                      ij                                                                        (                                                                          1                                      -                                                                              e                                                                                  -                                                                                                                                    Tw                                              p                                                                                                                                      RT                                                                                              2                                                ⁢                                                i                                                                                                                                                                                                                                                                                          )                                                                    ⁢                                                                      e                                                                                                                  -                                                                                  [                                                                                                                                    1                                                                                              T                                                                                                  2                                                  ⁢                                                  i                                                                                                                                                                                      +                                                                                                                                                                                            (                                                                                                                                                            D                                                      j                                                                                                        ·                                                                                                          G                                                      p                                                                                                        ·                                                    γ                                                    ·                                                    TE                                                                                                    )                                                                                                2                                                                                            12                                                                                                                                ]                                                                                                                    ⁢                                                                              t                                                                                  q                                          p                                                                                                                                                                                                                                                                          -                                                          g                                                              q                                p                                                                                                              )                                                2                                                                              )                                            +                                                          ⁢                                                          ⁢                                                                                                            W                      ⁡                                              (                        α                        )                                                              ·                                          m                      _                                                                                        2                                      }                          ,                            (        3        )            where p is the echo train index and P is the total number of different echo trains in a data acquisition sequence; wp is the weight applied the pth echo train. This minimization problem is usually solved by the least squares method, where m is the vector containing all the unknowns mij, i=1, 2, . . . I, j=1, 2, . . . J, and W(α) is a regularization matrix. One of the key procedures of the least squares method is to convert the expression (3) into the matrix equation (4)
                                                        ∑                              p                =                1                            P                        ⁢                          {                                                                    [                                                                                                                        (                                                          A                                                              L                                ×                                JI                                                                                      )                                                    p                          T                                                ·                                                                              (                                                          A                                                              L                                ×                                JI                                                                                      )                                                    p                                                                    +                                                                        W                          T                                                ⁢                        W                                                              ]                                    ⁢                                      m                    _                                                  -                                                                            (                                              A                                                  L                          ×                          JI                                                                    )                                        p                    T                                    ·                                                            g                      _                                        p                                                              }                                =          0                ,                            (        4        )            where L is the total number of echoes in the echo train; JI=J·I denotes the number of the unknowns to be inverted, and the vector gp contains the measured amplitude of the pth echo train. Due to the overwhelming large number of the matrix elements, the matrix operations performed in equation (4) is time-consuming, significantly slowing down the signal processing speed. Furthermore, processing using a huge amount of numerical operations, as is typical in such an analysis, can accumulate serious quantitative errors.
The second term in equation (3) is a regularization term. W(α) means that the entries of matrix W are functions of the regularization factor α (alpha). A simple example of W(α) is an identity matrix multiplied by α. The inclusion of a regularization term is important to many inversion problems using the least squares method, because least squares solutions are affected by data and rounding errors. The regularization dampens the error effects. It can be important to use adequate regularization such that it de-sensitizes noise influence while not introducing significant misfit. In general, an adjustable regularization coefficient is used. Normally, this adjustable parameter, α, is often determined by a cumbersome process involving multiple fitting trials with varying α.
Mathematical relation (3) can be expressed as,(A+W(α)TW(α))x=B  (5)
where x=[m11, m12, . . . , m1J, . . . , mI1, mI2, . . . , mIJ]T 
With the insertion of the regularization term, the solution vector x from equation (5) becomes a function of the factor α, denoted as x(α). Similarly, the error vector Ax−B also becomes a function of α, denoted as Ax(α)−B. If the norm of the error vector ∥Ax(α)−B∥ and the norm of the solution vector ∥x(α)∥ are used as the horizontal and vertical axis respectively to plot a curve, it will present an L-shape as shown in FIG. 1. The point on the curve with the maximum curvature corresponds to the optimal α. On the other hand, if a function of α is defined asd(α)=√{square root over (∥Ax(α)−B∥2+∥x(α)∥2)},  (6)then the plot of d(α) with respect to α can produce a parabolic-shaped curve, as shown in FIG. 2. The α value that minimizes d(α) corresponds to the stability of the inversion equation and the solution error. This α value is also referred to as the optimal α. Both the use of the L-curve and the use of the parabolic shaped curve can predict an optimal α based on the above criterion of the respective curves. Both methods may predict similar results. FIG. 3 shows an example where the maximum curvature of the L-curve (curve 342) and the minimum value of d(α) (curve 341) occur at almost the same α point. So as long as the correlation of ∥Ax(α)−B∥ with α and the correlation of ∥x(α)∥ with α are calculated using a list of trial α values, an optimal α can be determined either through the maximum curvature of the L-Curve or the minimum point of the parabolic type curve. Which of the above two methods, L-curve or parabolic shaped curve, is employed may depend on the performance of application in a given project. However, to obtain each point for the L-curve and for the parabolic curves requires performance of one iteration of the inversion process per point. As a result, calculation of multiple points in a curve takes a significant amount of time.